Trapping massless Dirac particles in a rotating saddle 
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We study particle motion in rotating saddle-shaped potentials. It is known that such rotating potentials can 
generate bounded motion for particles with a parabolic dispersion law through the balance of potential, centrifu- 
gal and Coriolis forces in the rotating frame. When applied to massless Dirac particles, for example electrons 
in graphene, such a potential is shown to lead to eigenstates that are spatially localized near the center of the 
saddle at certain energies. Although other states also exist at these energies, they are located further out from the 
center and this spatial separation tend to give the localized states a substantial life-time also when imperfections 
are present. 



It is well-known that it is impossible to trap a charged par- 
ticle with static external electric fields in three dimensions. 
Allowing for magnetic fields or time-dependent potentials 
changes this picture. Indeed, in atomic physics ions are rou- 
tinely trapped and studied in Penning [1,2] and Paul [3] traps. 
In a Penning trap the combination of an external perpendicu- 
lar magnetic field and a repulsive potential generate bounded 
motion in two dimensions (2D); when a particle tries to es- 
cape outwards it is deflected by the Lorentz force and pushed 
into bounded orbits. In the Paul trap a periodically driven 
quadrupole field is used to confine the ions in 2D. The physics 
of the Paul trap is often motivated by the closely related me- 
chanical problem of a point-like particle that is moving in two 
dimensions in a rotating saddle-shaped potential [3-5]. This 
problem (both the classical and the quantum version) can be 
solved exactly with elementary methods, as will be outlined 
below. In contrast the analysis of the Paul trap rests on the so- 
lutions of the Mathieu equation, which is considerably more 
involved. The relation between the rotating saddle potential 
and the Paul trap is akin to the relation between the Jaynes- 
Cummings model [6] and the Rabi model [7], in that the exact 
solution of the first is elementary while the analysis of the sec- 
ond one is considerably more involved [8]. As we shall see, 
the analysis of our system is greatly simplified by the fact that 
it is effectively time-independent in the rotating frame. 

A related and important problem in the field of graphene 
is the great difficulty of creating bound states by electrostatic 
means. The reason for this difficulty is that both electrons 
and holes have to be confined at the same time, but this is not 
possible in the absence of a (band) gap. Apart from the spe- 
cial case exactly at zero energy, where algebraically decaying 
states may be found [9, 10], essentially all states are extended 
scattering states, for a possible exception see [10]. Sticking 
to static electric fields, atomic-scale potentials are generally 
believed to be necessary for confinement. A fairly general 
classification of possible mass terms in graphene can be found 
in [11], examples that generate bulk gaps include static stag- 
gered sublattice potentials [12], antiferromagnetism [13], and 
spin-orbit coupling [14]. It is however presently not possi- 
ble to have sufficient control over these gaps. It is of course 
also possible to cut out a small flake and terminate the crystal. 
This creates confinement, but at the expense that the intrigu- 
ing Dirac fermion physics might be lost due to lattice effects. 



It would be better if smooth potentials could be used to create 
confinement. 

In this letter we present a scheme to trap electrons in 
graphene (and other systems described by a massless 2D Dirac 
equation) by subjecting them to a rotating saddle-shaped po- 
tential. This can be achieved by having two superimposed 
time-periodic quadrupole potentials with the proper rotation 
and phase shift. Our proposal is a natural extension of the re- 
cently popularized idea to modify material properties through 
spatially homogeneous periodic driving [15-19], to also allow 
for spatially dependent periodic driving. We present theoreti- 
cal arguments for trapped electron states by studying the Dirac 
equation in a rotating potential using semi-classical analysis, 
perturbative arguments, and a numerical investigation of the 
full quantum problem. To set the stage we first briefly review 
the behavior of the related quadratic problem, since this can 
be solved exactly and the result might be unfamiliar to non- 
experts. 

Classical particle. Let us begin by studying a classical par- 
ticle with mass m that moves in two dimensions in a saddle- 
shaped potential 

V{x) = ~(ulxl-u 2 2 xl), (1) 

that is rotating with angular frequency with respect to the 
lab frame, cji and u> 2 are frequencies that parametrize the 
shape of the potential. The classical equations of motion in 
the rotating frame are 

.ti = Q 2 xi + 2Qx 2 — UiXx, 

9 . 9 (2) 

&2 = Q X2 — 2ilil + U! 2 X2- 

On the right-hand side the first term is due to the centrifu- 
gal force and the second to the Coriolis force [20]. Solu- 
tions are easily found through the ansatz xi — cie~ tXt and 
x 2 = c 2 e~ lXt with ci and c 2 constants. A is determined from 
the solutions of the characteristic equation, which can be ex- 
pressed as A 2 — (fl±n 1 ) 2 -Q 2 l , w here fl 2 = (w 2 +w 2 ) /(4f2) 
and Sli = a/ (lo 2 — w 2 )/2 + tt 2 . For stability all A must be 
purely real, this gives two cases 

w 2 < uj 2 < ti 2 or 

ujI < u? 2 < Siof and uj\ < < ( ^ 9 + ^ . (3) 
1 " " 8(w 2 -w 2 ) 
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This means that when the saddle is more attractive than repul- 
sive the driving frequency has to be large enough to achieve 
stability, while if the saddle is more repulsive than attractive 
and not too repulsive there is a finite frequency interval where 
stable motion is possible. 

Quantum particle with parabolic dispersion. Let us now 
consider a quantum particle with a quadratic dispersion rela- 
tion that moves in two dimensions in the same rotating po- 
tential as in the classical case (1). In the rotating frame the 
stationary states are given by the solutions to the eigenvalue 
problem 



H p <P 



2m 



(4) 



Here L 3 = X\p 2 — x 2 p\ is the usual angular momentum op- 
erator. This Hamiltonian is straightforwardly diagonalized by 
a Bogoliubov transformation [21, 22], the result is 

H p = K\ + (a\a + - 1/2) - h\_{a ] _a_ - 1/2), (5) 



with \± — \/ ± f2i) 2 — fl|. Here a± are two indepen- 
dent bosonic annihilation operators. Stable motion is only 
possible when A± are both real, which gives the same sta- 
bility conditions as for the classical case (3). The complete 
set of square integrable eigenstates are the usual ladder of 
harmonic oscillator eigenstates built upon the vacuum annihi- 
lated by both a + and a_ . The fact that the energy eigenvalues 
are not bounded from below makes the system meta-stable in 
the presence of perturbations, the life-time can nevertheless 
be quite long [4] . 

Semi-classical treatment of massless Dirac particles. Af- 
ter this introduction we now consider the 2D massless Dirac 
Hamiltionian in an external potential 

ih^lab = VPlab ■ CT^lab + V(Xj o6) t)^i ab , (6) 

here v > is the velocity of the particles. In the graphene 
problem this is the Hamiltonian close to one of the K-points 
and the Pauli matrices <x act in sublattice space [23]. Sta- 
tionary states can be found by going to the rotating frame 

^lab{t,Xi ab ) = e -iEt/h e ~i(nt/2)a 3qi ^^ which j ea( j s tQ ^ 

eigenvalue problem 



H^> = vp ■ cr* + y(x)* - OJ 3 * = E^f. 



(7) 



Here J3 = L3 + ha^ /2 is the total angular momentum opera- 
tor, which has eigenvalues that are half-integer multiples of fv. 
±fi/2, ±3fi/2, ±5fi/2, etc. Since there is no mass in the bare 
graphene problem we define a dynamically generated mass 
via mi) 2 = fi|fi|/2 and continue to use the parametrization 
of the potential in (1). Let us now perform a semi-classical 
analysis of this Hamiltonian. Diagonalizing the matrix struc- 
ture we get two semi-classical energy bands corresponding to 
electrons (+) and holes (-): E± = ±E p - Qi 3 + V(x), 
with E p = \J v 2 p 2 + (mn 2 ) 2 . The semi-classical equations 



of motion for wave-packets in multi-band systems such as this 
one are modified by a Berry curvature term [24, 25] 



x = X7 P E± - p x Sl±(p)/H, 
p = -V x £±. 



(8) 



In our case the Berry curvature is Q± — ±h 3 v 2 H,n 3 /(4E 3 ). 
Eliminating x from the equations of motion we arrive at 



Pi 



n{u;j-uj 2 ). n 2 co 2 h\n\ujl h 3 n\n\uj 2 . 



-P2 



W 2 " PlT 2EL 



Pl± 



8E p 



-P2, 



P2 = 



^(w|-w 2 ). ft 2 w| ,h\n\u)% n 3 n\n 



-pi 



P2 ± 



2EL 



-P2 ± 



8E p 



-Pu 



(9) 



which determines the semi-classical wave-packet dynamics 
for electrons (upper signs) and holes (lower signs). For |p| <C 
mv we can linearize these equations and look for stable solu- 
tions exactly as for the classical particle above, i.e., we make 
the ansatz pi = cie~ lXt and p 2 = c-2,e~ %xt and solve for A. 
Note that the linearization implies that E p s» fi|fi|/2 = mv 2 . 
Demanding that all of the eigenvalues are real for wave- 
packets constructed in both of the bands we get the necessary 
condition (V2 — 1/2)cj 2 < lo 2 , < w 2 where w 2 (w 2 ) is the 
greater (smaller) of w 2 and ui 2 . As a result the two curva- 
tures {J 2 and uj 2 are not allowed to be too different. In addi- 
tion there is the constraint w 2 < f2 2 < uj 2 > f c (\u! < /uj > |), with 
f c ( x ) = 2 ; 2 /{(l-a ; 2 ) 1 / 3 [(l + x) 2 / 3 -(l-2;) 2 / 3 ] 2 }. This im- 
plies that if w 2 ^ (J 2 the motion is only stable for a finite inter- 
val of driving frequencies 57. Specializing to the particle-hole 
symmetric case (see discussion below) where uj\ = cj| = uj 2 
the solutions are A 2 = Q 2 and A 2 = Q 2 — lj 4 /Q 2 , and 
hence the sufficient condition for stability in this case is sim- 
ply n 2 > uj 2 . Although this analysis hints at what to expect 
from a full quantum mechanical solution it should be viewed 
with some skepticism. For the semi-classical analysis to be 
valid the time scale associated with the wave-packet dynam- 
ics should be much longer than the time-scale of the gap, i.e., 
|A| <C This is never true for all of the semi-classical 
bounded solutions, although it is for some. 

Performing the analogous analysis for an ellipsoidal poten- 
tial (taking uj 2 < and ui\ > 0) we find stable semi-classical 
solutions for sufficiently small driving frequencies. In this 
case, however, we always find that |A| > |fi| so that these 
solutions can never be trusted. 

Quantum treatment of particles with linear dispersion. To 
study the quantum problem it is convenient to first reformulate 
it in dimensionless variables. Denoting uj 2 = (u 2 +u> 2 )/2 and 
A 2 = uj 2 — u 2 we introduce the characteristic length scale 
I = 2v/u and write Xj/l = Xj. We also go to a different 
gauge ^ — > exp(isQX 1X2)^ with sq = sign(17) which takes 
Pi — > Pi + sq$2 and i?2 ^ > P2 + snxi so that the Hamiltonian 
goes into 



H 

m 



-j 3 +^(p. ff +« SJ i-<TV 3 )+^(i 2 +i 2 ). (10) 
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We now reexpress H in an oscillator basis where Xj and dj 
are represented with creation and annihilation operators in the 
usual way — (dj + Xj)/y/2. In terms of the bosons b l = 
(—ia 1 + a 2 )/v2 and b 2 = (— ia, — a 2 )/v / 2 we get 



H 
Ml 



U) 



o 

h.c. 







2VL 



snA 2 _ 2 

r 

2uj 2 

b l - 
h.c. 



with 



- J 3 = b\by - b\b 2 - <7 3 /2, 



r 2 



r,2 



~2 



= bib, 



blb 2 + 1-6^2 



(11) 



(12) 



When A 2 = the spectrum is particle-hole symmetric. This 
can be implemented by the transformation b\ — > i& 2 and 
b 2 — » — «6i combined with the rotation i? = e~ 47rcr3 / 4 <7 2 on 
the Hamiltonian (11), this also takes J 3 — > — J 3 . 

We note that all of the terms on the first line of (11) 
commutes with J$ which has half-integer eigenvalues j — 
±1/2, ±3/2, . . ., whereas the term on the second line cou- 
ple states with values j to j ± 2. If the first term is the largest 
one a straightforward degenerate perturbation theory calcula- 
tion to second order in the remaining terms, similar to the one 
in [16], gives an effective Hamiltonian that commutes with J 3 



Hefi 

~Ml 



20 



per 



/ LO \ 2 f 2 



<?3 



snA 2 . 2 

r 

2cj 2 



(13) 



This is a 2D Dirac equation with a spatially varying mass term, 
which is an established way to confine Dirac particles [26]. 
From H e ff we also get a simple necessary condition for con- 
finement by demanding that the strength of the mass term is 
larger than that of the parabolic potential. This implies that 
il 2 |A 2 | < w 4 /4 and puts an upper bound on il 2 that goes 
to infinity in the limit A 2 — > 0. On the other hand we also 
need lu <C \il\ for the perturbative treatment to be valid. This 
result is therefore consistent with the semi-classical stability 
analysis above with an upper and lower bound on the allowed 
driving frequency, but less numerically precise. Solving H e g 
numerically for j = —1/2 we find that the perturbative treat- 
ment is self-consistent for states with energies close to Ml/ 2. 
The results obtained from _ff e ff are also in approximate agree- 
ment with those of the full Hamiltonian H. 

For simplicity we now consider the particle-hole symmetric 
case A 2 = and that fl ^S> lu. When to = the eigenener- 
gies are —h£lj with j = ±1/2, ±3/2, . . ., all of which are 
infinitely degenerate. Perturbatively the degeneracy is split 
when lu ^ 0. It is straightforward to write down the matrix 
representing the Hamiltonian (1 1) in the standard occupation 
number basis on each sub-lattice: {\ni, n 2 )} Illi?l2 >o defined 
by b\b 1 \m,n2) = ni\ni,n 2 ) and b\b 2 \ni,n 2 ) = 7i 2 |rti,n 2 ). 
Truncating the matrix with n 1: n 2 < N c and diagonalizing 
it numerically we find that some eigenstates and the associ- 
ated energies converge near E « Ml/ 2, see Fig. 1. These are 
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FIG. 1. Normalized energy eigenvalues of H/mv 2 closest to 1 ob- 
tained by numerical diagonalization of the full particle-hole symmet- 
ric quantum problem for uj/Q. — 0.1 and different cut-offs JV C . The 
corresponding matrix sizes are Some eigenvalues, cor- 

responding to j = —1/2, converge and remain stable as N c is in- 
creased. For sufficiently large N c (~ 100 in this case) new states 
with energies near E/mv 2 — 1 appear and that are also stable upon 
further increase of the cut-off. These have a large overlap with states 
of angular momentum — 1/2 ± 1. 



E= 1.100 mv 2 M ||l III 




E= 1.088 mv 2 




E= 1.037 mv 2 


E= 1.001 mv 2 1 
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FIG. 2. Weights in different oscillator states for some of the eigen- 
functions corresponding to Fig. 1 (N c = 250). Weights for different 
energies are shifted and rescaled independently. Blue and red (in- 
verted for clarity) denote different sublattices. The states with weight 
close to the origin are to ~ 99% built from j — —1/2 states and the 
other states have approximately equal weight inj = — 1/2 ± 1. 



the localized states of interest, and they have a large weight 
in the j = —1/2 sector. There are also other states near this 
energy, but they are constructed from large values of t%\ and 
7i2. Increasing the cutoff N c the first such states that appear 
have large weight in the j = — 1/2 ± 1 sectors. A simple 
estimate of the associated quantum numbers is obtained by 
the approximation b\ ~ ^fn\ ~ b 2 ~ ^/n^ and therefore 
ni ~ n 2 ~ (il/oj) 2 . Some of the numerically obtained wave- 
functions are shown in Fig. 2 and are consistent with this es- 
timate. Clearly more states appear when larger cut-offs are 
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FIG. 3. Sketch of an experimental setup (top view) that can be used 
to generate the rotating saddle potential. Two sets of mechanically 
fixed electrodes (a and b) are placed symmetrically around the center 
of the sample. Time-dependent potentials Q are applied to the elec- 
trodes with the indicated strengths and relative signs so that each set 
generate a quadrupole field in the center. When the amplitudes be- 
tween the to sets have a 7r/2-phase shift, e.g. Q a (t) ~ Qcos(2S!t) 
and Qt(t) = Q sin(2f2£), the potential in the center V(x; a j,, t) be- 
comes that of a rotating saddle. 



considered corresponding to even bigger \j | . An important ob- 
servation is that since the different sets of wavefunctions have 
essentially non-overlapping support, smooth disorder poten- 
tials will only weakly couple the different sets, which implies 
that the life-time of the localized states is large. 

Possible application in graphene and on surface states of 
3D topological insulators. A rotating saddle potential is in 
principle easily be generated by a superposition of two oscil- 
lating quadrupole fields, as sketched in Fig. 3. It is interesting 
to note that, in the rotating frame, a localized disorder po- 
tential on average generate a rotationally invariant potential 
profile. If the driving frequency is large enough this first or- 
der result is the dominating contribution. We therefore expect 
that the bound states are only weakly affected by disorder and 
imperfections. 

Experimentally the localized states can be probed using 
scanning tunneling spectroscopy or via transport through the 
bound states. For clear experimental signatures of the trapped 
states the typical energy separation between neighboring j- 
values near the saddle center H\Q\ should at least exceed 
the temperature scale T. Experimentally frequencies above 
~ 10 GHz are very difficult to achieve with conventional ra- 
dio frequency methods, limiting the temperature range of such 
setups to the sub-Kelvin scale. In the optical frequency range 
the temperature is no longer an issue. In this case the design 
of an appropriate laser setup to achieve the desired inhomoge- 
neous field configuration is however more challenging as the 
dipole approximation is not accurate enough. 



To summarize we propose to use a rotating saddle-shaped 
potential to trap massless Dirac particles in 2D. This has po- 
tential applications in graphene, three-dimensional topologi- 
cal insulator surface states, and cold atom systems. Our study 
also demonstrates the richness of particle dynamics in external 
time- and space-dependent fields. 
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